(N 

O 
(N 



Three particles in a finite volume: 
The breakdown of spherical symmetry 

Simon KreuzeiQ and Harald W. GrieBhammeilll 

Institute for Nuclear Studies, Department of Physics, 
The George Washington University, Washington DC 20052, USA 

(Dated: March 3, 2013) 

Abstract 

Lattice simulations of light nuclei necessarily take place in finite volumes, thus affecting their 
infrared properties. These effects can be addressed in a model-independent manner using Effective 
Field Theories. We study the model case of three identical bosons (mass m) with resonant two-body 
■ interactions in a cubic box with periodic boundary conditions, which can also be generalized to the 

three-nucleon system in a straightforward manner. Our results allow for the removal of finite volume 
effects from lattice results as well as the determination of infinite volume scattering parameters from 
the volume dependence of the spectrum. We study the volume dependence of several states below 
the break-up threshold, spanning one order of magnitude in the binding energy in the infinite 
volume, for box side lengths L between the two-body scattering length a and L = 0.25a. For 
CJ . example, a state with a three-body energy of — 3/(ma^) in the infinite volume has been shifted to 

^ i — 10/(ma^) at L = a. Special emphasis is put on the consequences of the breakdown of spherical 

symmetry and several ways to perturbatively treat the ensuing partial wave admixtures. We find 
their contributions to be on the sub-percent level compared to the strong volume dependence of 
• the S-wave component. For shallow bound states, we find a transition to boson-diboson scattering 

. behavior when decreasing the size of the finite volume. 
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I. INTRODUCTION 



Quantum chromodynamics (QCD) is the theory underlying strong interactions. However, 
ab initio calculations of hadronic and nuclear properties remain one of the largest theoretical 
challenges of the Standard Model. Lattice simulations provide a numerical approach, but do 
at present not usually operate at the physical point (for reviews on the lattice simulations 
of light nuclei see, e.g., [H, B])- In a complementary approach. Effective Field Theories 
(EFT) describe the effective degrees of nuclear physics, namely nucleons and pions, and 
allow for accurate calculations of low-energy observables with a direct link to QCD through 
the symmetries of the theory, see e.g. 

In lattice simulations, the QCD path integral is evaluated in a discretized Euclidean space- 
time using Monte-Carlo simulations. This approach requires a large numerical effort, which 
in turn strongly constrains the parameters of the simulation. In particular, the considered 
system is necessarily placed in a finite volume. Present day calculations use cubic boxes 
with periodic boundary conditions and relatively small side lengths of a few fm. Momentum 
quantization due to the boundary conditions causes a shift of the finite-volume spectrum 
relative to the infinite- volume energies |9|] . A model-independent determination of this shift 
is necessary in order to extract physical observables from lattice results. Calculations at large 
volumes are possible, but rendered inefficient by the enormous numerical effort necessary. 
Therefore, the use of well-known physics to perform the extrapolation is warranted. 

The volume dependence of the spectrum also provides access to scattering parameters. 
Most prominently, Liischer showed that infinite volume scattering phase shifts as well as 



resonance properties are encoded in the finite- volume spectrum of two-particle states [10|, lU 



The correlation function for the three-nucleon system in the triton channel has been 
calculated in Lattice QCD recently [l2| , but because of the relatively large uncertainties no 
triton properties could be extracted. With quantitative lattice data on light nuclei and their 
scattering properties within reach, Liischer's results need to be extended to the three-body 
sector in order to understand the finite volume effects in these results. 

The desired incorporation of well-known physics, combined with the demand of model 
independence and the ability to achieve a given level of accuracy, are preconditions well 
met by the EFT approach. In the nuclear sector, the pionless EFT, which is valid for 
processes with typical momenta below the pion mass, has been successfully used to describe 
the properties of light nuclei (see, e.g., [1, 0, 0] for reviews). To leading order, the three-body 
sector is described by the nucleon-nucleon scatterin g le ngths in the ^5*0 and channels 



and a Wigner SU(4) symmetric three-body force [l3Hl5 



The pionless EFT in the nucleonic sector can be seen as part of a larger family of EFTs for 
systems with resonant two-body interactions. They are characterized by the appearance of a 
two-body scattering length a that is unnaturally large compared to the range of interaction. 
Such systems display interesting universal properties. If a is positive, two particles of mass m 
form a shallow two-body bound state with binding energy l/(ma^), independent of the 
detailed mechanism generating the large scattering length. For example, ^He atoms have 
the unnaturally large scattering length 189 oq, where oq is the Bohr radius. Using the 
aforementioned relation to extract a from the energy of the helium dimer, on the other 
hand, yields a scattering length of 182 ao, deviating by only 3.7%. In the pionless EFT, the 
shallow two-body bound state is identified as the deuteron. 

In the three-body system, the universal properties include the Efimov effect [l5|. If at 
least two of the three particles have a scattering length \a\ that is large compared to the 
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range of their interaction, a sequence of three-body bound states exists. In the hmit of 
diverging scattering length, there are infinitely many geometrically spaced bound states 
with an accumulation point at threshold. This spectrum is the signature of the Efimov 
effect, namely a discrete scaling symmetry whose consequences can be calculated in an EFT 
for short-range interactions. Here, it appears because the renormalization group flow of the 
three-body coupling is a limit cycle (see, e.g., [l^, 17| for a review). This is the case for the 



Wigner SU(4) symmetric three-body force in pionless EFT [l3j, but also for the three-body 
force needed to renormalize the simpler EFT of three identical bosons 18|]. In this work, 
we will therefore study systems of three identical bosons inside a cubic box with periodic 
boundary conditions using the EFT framework. Transferring the methods developed in this 
publication to the pionless EFT will give an understanding of the finite volume corrections 
in lattice calculations of the three-nucleon system. 

In this paper, we study the volume dependence of three-boson states below the three-body 
break-up threshold inside a cubic box with periodic boundary conditions. One of the authors 
performed similar studies before We developed the new framework presented in this 



paper for several reasons. First, the new approach provides access to the energy region above 
the boson-diboson break-up threshold. Second, it does not rely on the use of basis functions 
but remains in close contact with the discrete space of allowed momenta. It further allows for 
a faster numerical implementation when higher partial waves are taken into account. Thus, 
for the first time, we are able to study the size of these higher partial wave contributions in 
a systematic manner as well as how they can be included perturbatively. 

The application of EFT to finite volumes, especially for three-particle systems, has 
sparked interest in recent years. The volume dependence for three spin-1/2 fermions with 
perturbative interactions in a box has been studied previously ^3|. Epelbaum and collabora- 
tors have calculated the energy of the triton in a finite volume by implementing a discretized 



version of chiral EFT on a lattice [2J]. Also, the triton has been considered in pionless EFT 



in a nuclear lattice formalism but the volume dependence was not investigated [25|. Other 



volume shapes have also been studied, most prominently the harmonic oscillator [26|, |27 
which allows to compute finite volume effects inside atomic traps. 

This paper is organized as follows. In the next Section, we derive the basic equations 
and sketch the ideas of our numerical implementation. We furthermore present two different 
approaches to include higher partial wave corrections perturbatively. In Section IIIIl we 
provide a detailed discussion of our numerical results, including a study on the perturbative 
nature of higher partial waves. Section HVl summarizes our findings and provides an outlook 
for future work. 



II. FRAMEWORK 



In the following, we derive a set of coupled equations governing the partial waves of the 
amplitude and explain how to prepare these equations for a numerical implementation. In 
the last part of this section, we show two ways to take higher partial waves into account 
perturbatively. 
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FIG. 1. Integral equation for the boson-diboson amplitude, represented by the shaded circle. The 
single line denotes the single boson, while the double line denotes the diboson field. 



A. Sum equation for the amplitude 

The Lagrangian for three identical bosons interacting via short-range forces can be written 



as (see, e.g., |l7l. Il8j) 



£ = (tdt + ^ + ^dU - ^ (dV + H.c.) - + • • ■ , (1) 

\ 2 / 4 4 6b 

where the dots indicate higher order terms of the effective theory. The Lagrangian is formu- 
lated in terms of the boson field ip and a non-dynamical auxiliary field d with the quantum 
numbers of two bosons. Units have been chosen such that h = m = 1, where m is the mass 
of a single boson. 

The system of three bosons is assumed to be contained in a cubic box with side length L 
and periodic boundary conditions. This leads to quantized momenta ^n, n & 1? . As a 
consequence, the loop integrations of the infinite volume case are replaced by sums. If loop 
sums are divergent, they are regulated by a cutoff A similar to the infinite volume case. 

While the finite volume modifies the infrared regime of the theory by introducing the 
low- momentum scale 27r/L, the ultraviolet behavior does not change. Therefore, its renor- 
malization is the same in the finite and infinite volume cases. Of course, this statement is 
only valid as long as the infrared and ultraviolet regime of the theory characterized by the 
scales 2tt /L and A, respectively, are well separated, i.e., AL ^ 1. We will explicitly demon- 
strate that the numerical results in this publication are indeed renormalized. The cutoff 
dependence in the two-body sector can be removed completely by matching the coupling 
constant g2 to a low-energy two-body observable, namely the two-body scattering length a 
or, if applicable, the two-body binding energy. As the cutoff is increased, the three-body 
coupling approaches a renormalization group limit cycle, leading to a log-periodic depen- 
dence of the coupling constant on the cutoff A. For convenience, is often expressed in 
terms of a dimensionless function H via g^ = — 9(?| if(A)/A^. The phase of the log-periodic 
dependence of H has to be fixed from a three-body datum (l6l-[l8l.l28j. 

The central quantity in the three-body sector is the boson-diboson amplitude. The deriva- 
tion of the equati on g overning this object in finite volume has been shown in detail in a pre- 
vious publication [20, 2^ and is therefore only briefly sketched in the following. A discussion 



of the infinite- volume case can be found, e.g., in [17|, |28j. Starting point for both cases is the 
Lippmann-Schwinger equation for the boson-diboson amplitude. In a finite volume, there 
is no continuous scattering spectrum, but only discrete energy levels. The analytic struc- 
ture of the amplitude is given by a number of simple poles at the energies of these discrete 
states, just like the analytic structure of bound states in the infinite volume. In addition, 
the incoming and outgoing quantum numbers separate. The residue at the pole is the ampli- 
tude J-". Comparing the residues of both sides of the Lippmann-Schwinger equation at any 



4 



of the discrete energy levels yields the homogeneous equation depicted diagrammatically in 
Figure [H 

In the center of mass frame, the momenta of the outgoing single boson and diboson are p 
and — p, respectively. The outgoing single boson is taken to be on-shell, giving the four- 
momentum (p^/2,p). The four-momentum of the outgoing diboson line is [E — — p), 
making the total kinetic energy E a parameter of the equation. The equation depicted in 
Fig. [T]then serves as a consistency condition. Values of the energy parameter E for which the 
equation has a solution are identified as the energies of the discrete finite volume spectrum. 

Using the kinematics outlined above, the sum equation for the boson-diboson amplitude J-" 
reads 

^^P) = % E Z{E:P.q)r{E-q)m- (2) 



The quantity Z{E;p,q) contains the interaction kernel of the sum equation, while T{E]q) 
is essentially the diboson propagator inside the loop. 

For the diboson lines in the three-body equation depicted in Fig.[Tl the interacting diboson 
propagator has to be used. This quantityjs obtained by dressing the bare propagator, which 
is a constant, with bosonic loops (cf. (l8j). In the sum equation ([2]), it appears via the 
quantity pj, i22)] 



.(.;,)^ i-V^^-.^Ezj^e-l.^^-, (3) 

which is identical to the full propagator up to a pref actor. The absolute value of the integer 
three- vector J is denoted In the limit L — )■ oo, the propagator reduces to the infinite 
volume expression. 

The interaction kernel is given by the one-boson exchange term and the boson-diboson 
contact interaction. Since the time direction is considered to be of infinite extent, the 
integration over the loop energy remains continuous and can be performed analytically by 
using the residue theorem. The resulting interaction kernel depends only on incoming and 
outgoing three-momenta as well as the total energy: 



Z{E-p,q) 



1 ^ H{A) 



p^+p-q + q^ — E 



(4) 



In the present calculation, only energies below the threshold for the break-up into three 
individual bosons, i.e. < 0, are considered. This keeps the investigated states below any 
inelastic thresholds, in particular the Nd — )■ NNN threshold. 



B. Methodology 

The amplitude J-" is strictly speaking only defined for the discrete momenta inside the 
box. One ansatz for the solution of the sum equation ([2]) is therefore to use the quantized 
momenta themselves as sampling points. In this approach, the value of the amplitude at the 
sampling points is written as fi = J^{pi), where the index i labels all momenta pi G 
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TABLE I. Coefficients CAiim of the kubic harmonics [3l[ | as weh as -STai^ in Cartesian coordi- 
nates [2a] for £ = 0, 4. 



with \pi\ < A. This yields the finite-dimensional eigenvalue equation 



h (5) 



The range of box side lengths that allow for the use of this ansatz is limited. The number 
of points to consider roughly scales as (AL)^, while the complexity of the necessary matrix 
diagonalization scales like the number of points to the third power. We will therefore only 
show results for L < a in this publication. Another drawback of this approach is that it is 
cumbersome to disentangle the contribution of different partial waves to the finite volume 
spectrum. We will refer to this ansatz as the Grid approach. 

The momentum quantization inside a finite volume is tantamount to a reduction of the 
spherical symmetry of the infinite volume to a discrete symmetry, namely in our case the 
point symmetry of the cube. In the language of group theory, the infinitely many irre- 
ducible representations of the rotational group SO (3) become reducible in terms of the five 
irreducible representations of the cubic group O. Hence, a quantity ips transforming accord- 
ing to the irreducible representation s of O can be written in terms of the basis functions of 
spherical symmetry, i.e. spherical harmonics Y^m, via 

Up) = Y, RMKsitiP), with K^etiP) = J2 CtLyUP). (6) 

where p = p/p = {piip2iPz) is the unit vector in j>-direction and Ru is the radial function. 
The index t is needed if the representation labeled by £ appears in the irreducible represen- 
tation s more than once. The coefficients Cf^^ are real and normalized for given values of 

— 1. The linear combinations Ksu of spherical harmonics are 



s, £ and t via 



called "kubic harmonics" fsic!) l29i. The values of the coefficients Cf^^ are known for i as 



large as 12 [30[ and are readily computed from group theory [31|. In this work, we restrict 
ourselves to amplitudes J-" transforming under the trivial {Ai) representation of the cubic 
group, because it contains the trivial representation of spherical symmetry, namely i = 0. 
This is because the investigated Efimov states are S-wave states in the infinite volume. In 
addition, Ai receives contributions from £ = 4, 6, 8, . . . Since all ^-values appear only once 
below i = 12, the index t is dropped in the following. The coefficients C^i^m for £ = 0,4 are 
summarized in Table [T] together with the resulting kubic harmonics in Cartesian coordinates, 
normalized to An. 

The expansion of the angular dependence of the amplitude in spherical harmonics allows 
for an assessment of the impact of higher partial waves or, in other words, of the loss of 
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spherical symmetry in the finite volume. Moreover, it allows for an analytic calculation 
of the angular dependence, leaving only the radial part for a numerical treatment. This 
significantly reduces the complexity of the numerical problem at hand compared to the 
Grid approach. We therefore determine the p-dependence of the right hand side of Eq. ([2]). 
The only quantity where p is present is the one-boson exchange part of the interaction 
kernel Z{E]p,q), given by the first term in Eq. (jlj). Consider therefore the quantity 



f d^p 
'52 An 



KAit 



.Ceo 



1 



1 



a + p-Cq a — p-Cq 



(7) 



where we have set a = p^ + — E for simplicity. The summation of all elements of the 
cubic group O makes the expression in square brackets invariant under any cubic rotation 
of p and q. It therefore transforms according to the trivial representation, denoted Ai, of the 
cubic group. The inclusion of the second term makes it moreover invariant under additional 
parity transformations. The expansion coefficients of this quantity in kubic harmonics are 
given by Ip\q). The analytic evaluation of the angular integral is shown in Appendix Rl 
The result is 

I^'\Q) = -Qi(-)KA.iiQ), (8) 
P pq \pqj 

where Qp is the £th Legendre function of the second kind as defined in j32|. Special care has 
to be taken of the limiting cases p, g — )■ 0. The result /o(^ = -^p(O) = coincides with 
the zero- momentum limit of ([8]). 

We will now use this result to simplify - at least from a numerical viewpoint - the sum 
equation (|2]). The summation over the discrete set of three- vectors can be rewritten in terms 
of a summation over the elements of the cubic group as follows: 

q-e^zs '?6^(Z3) Ceo 

where {1?) = G : gi > 52 > 93 > 0}. All other vectors in 1? can be generated by 
applying cubic rotations and the parity operator to the vectors in (Z^), i.e.. 



U (C + CP) = Z\ (10) 

Ceo 

However, the summation over all cubic rotations in ([9]) might introduce double counting. 
Specifically, if the original vector contains identical entries or zeroes, several cubic rotations 
have the same image vector. This is accounted for by the factor mul(g)/48, where mul(g) 
is the number of vectors that contain the same entries as q up to ordering and signs. The 
identity Eq. ([9]) was also used to reduce the dimensionality in the numerical implementation 
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of the Grid approach. Applying ([9]) to Eq. ([2]) yields 



(^i) 



L3 



C60 



Stt 

Z3 



E 



1 



(^l) 



9 9 



48 
-E 



X mu.l{q)T{E; q)J^{q) 



KAAp)KAAq) + ^^KaA 



(11) 



where the £-sums run over the partial waves contained in the Ai-representation, denoted by 
the upper summation limit (v4i), and J^{Cq) = J^{q) E O was used in the first line. 

Comparing the coefficients of Ka_Ap) Eq. ffTTj) yields a set of coupled equations for 
the radial functions Re. Let Ru = Re{pi), where the pi are the possible absolute values 
of vectors in smaller than the cutoff A. The resulting expression has the form of an 
eigenvalue problem 

Ru = Mei.e/j{L; E) Rpj, (12) 

where the matrix M is given by 



M, 



87r \ - 



iPiQj V Pi(lj 



E\ H(A) 



A2 



(13) 



X mul(g)r {E; qj) KA^e'{q). 



We use a root finding algorithm to determine values of the energy parameter E such that 
the eigensystem (fT2!) contains an eigenvalue 1, yielding the discrete energy levels inside the 
finite volume. 

For the discussion of the perturbative approaches, it will be convenient to write the matrix 
equation ( !T2|) in block matrix form in angular momentum space, viz. 







1 Moo 


Mo4 






M40 


M44 













Ra 



(14) 



with {Re)i = Ru and {Mui)ij = M^i-ej. The dependence of the block matrices on the box 
side length and the energy parameter is suppressed in this notation. The off-diagonal block 
matrices in are responsible for the mixing of different partial waves. 

We will refer to the method described above as the Sum approach. The conceptual advan- 
tage of this approach is the possibility to disentangle the contributions of the different partial 
waves to the energy shift in the finite volume. On the numerical side, the dimensionality 
of the eigenvalue problem is greatly reduced in comparison to the Grid approach: if there 
are several three-momenta with identical absolute value, they all increase the dimensionality 
of the Sum approach only by one. The differences in dimensionality and runtime will be 
discussed at the end of Section IIIII 
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In earlier works 19, 20, 22l| the loop sum was rewritten into a sum of integrals via 
the Poisson equation. This explicitly recovers the infinite volume form of the amplitude. 
This is appropriate for bound states, where this approach was used. It is, however, not 
appropriate when extending the formalism into the energy region of elastic boson-diboson 
scattering. Here, the infinite and finite volume amplitudes have fundamentally different 
analytic structures, namely a scattering continuum and a series of poles corresponding to a 
discrete spectrum, respectively. This energy region is accessible with the new Sum approach. 
We will refer to the older approach as the Poisson approach. 



C. Perturbative approaches 

In the infinite volume, bound states and low-lying scattering states are predominantly S- 
wave states. In the finite volume, there are admixtures from higher partial waves stemming 
from the breakdown of the spherical symmetry, as discussed above. For reasonably large 
volumes, these admixtures are expected to be small compared to the S-wave only part. 

To see that this is indeed the case, consider the argument of the Legendre function Qi in 
the definition of the matrix elements (IT^ . The argument {p'f + q"^ — E)/ (piqj) is guaranteed 
to be larger than 1 (remember that E < here) and will be large if the ratio Pi/qj is either 
very large or very small. For large arguments z ^ 1, the Legendre functions Qe{z) scale as 
z~^~^ (see, e.g.,j32|). Therefore, 



Qo{z) 



"^ if z > 1, (15) 



which shows the suppression of the G-wave contributions. The momenta Pi, qj are multiples 
of the low-momentum scale 211 /L. Consider, for example, the pair of momenta with the 
largest possible ratio, namely A and With these choices, and the condition AL 3> 1 from 
renormalizability, the argument of the Legendre function becomes ^(A -|- (£^3 /A)). Using 
the expansion of the Legendre functions for large arguments, we can estimate the scaling of 
the G-wave correction by 

The scaling of other matrix elements with the box side length is not as drastic but still 
governed by the strong suppression evidenced in Eq. ( |T5l) . 

A second suppression mechanism comes from the kubic harmonics in the definition of the 
matrix elements (fT3|) . While the kubic harmonic Ka^o is just a constant, the corresponding 
functions for £ 7^ show sign oscillations (cf. Table [T]). These lead to cancellations among 
the terms from vectors with the same absolute value. These corrections are themselves inde- 
pendent of the size of the volume as they refiect the breakdown of the spherical symmetry. 
However, there are more cancellations among integer vectors of high absolute values, where 
the number of vectors with identical absolute value is large. These vectors are reached for 
AL ^ 1, leading to more cancellations within the summation over the possible momenta. 

Overall, the two suppression mechanisms in the definition of the matrix elements intro- 
duce the generic hierarchy 

Moo > Mo4 > M40 > M44 > Meo > . . . (17) 
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In the following, we discuss two different approaches to take advantage of this ordering. The 
first one can be described as a partially resummed calculation, while the second approach 
is developed as a strict perturbative expansion. In the following discussion, we will use the 
notation introduced in Eq. (III]) . 

For the Partial Resummation approach, consider the second row of Eq. (1141) . 



Ra 



MaoRo + MaaRa + 



(18) 



Inserting this expression for i?4 into the equation for Rq given by the first row of Eq. ([T] 
yields the term Mo^M^^R^ containing an additional factor of KAi^iq) compared to M04M40R0. 
This term is therefore discarded as parametrically small in the partial resummation ap- 
proach. More generally, the partial wave amplitudes Ri for £ 7^ are approximated by 
Ri^o = MiqRq, or equivalently by setting all M^^q/vo zero: 



Ra 



( 



Moo 


Mo4 




/ 


Rq 


M40 








R4 










/ 


\ • 



Roc^Imoo + Y, ^o^^^o ) ^0, (19) 



where the sum over i is truncated at £ = or £ = 4 in the calculations presented in this 
paper. This is possible due to the general hierarchy f lT7|) . Using this equation, the energy 
of the state can be obtained by the eigenvalue method described above. 

For the Strictly Perturbative approach, we first establish how to obtain the perturbed 
energies in the eigenvalue method. Consider the general eigenvalue equation 



T{E) = K{E)T{E) 



(20) 



for an amplitude T and a matrix K, where both quantities depend on a parameter E that is 
tuned such that T is an eigenvector of K with eigenvalue 1. If the matrix K contains a small 
perturbation K[E) = Kq{E) + eKi{E) with e <^ 1, both T{E) and E can be expanded in 
powers of e, retaining only the linear term: 



T{E) = To{E) + eTi{E) + 0{e^) 
E = Eo + £Ei + 0{e^), 



(21) 



such that 



To{Eo) = Ko{Eo)n{Eo) (22) 
holds. Furthermore, the ii^- dependence of Kq and To can be linearized around Eq as 



Ko{E) = Kq{Eo) + eE,K',{Eo) + 0{e^) 
To{E) = To(Eo) + eE,T;,{Eo) + 0{e'), 



(23) 



where the prime indicates differentiation with respect to the energy parameter. Inserting 
the expansions ( 12T|) and (!23|) into (!20|) and comparing orders of e yields to order just 
Eq. (122 p . To order e^, we obtain 



EiT^ + Ti = E,K'oTo + K,To + Ko[E,T;, + T,] 



(24) 
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where all quantities are evaluated at the leading order energy Eq. Multiplying by Tq{Eq) 
from the left hand side and using the daggered version of Eq. fl22|) yields the desired equation 
for the energy shift in lowest order in e 



E, = -Zo^. (25) 



Note that the leading order amplitude Tq appears only at the leading order energy. 

In the present framework, the leading order matrix Kq is given by the block matrix Mqq 
in Eq. (HM . while the leading order amplitude is accordingly the amphtude Rq determined 
such that 

Ro{Eo) = Moo{Eo)RoiEo) (26) 

holds, cf. Eq. (1221) . Since the dependence of Mqo on the energy parameter is known from 
Eq. (fT3|l . the differentiation in the denominator of (125|) can be carried out analytically. As 
can be read off from Eq. ( !T9|) in the discussion of the partial resummation approach, the 
perturbation is given by 

Ki = J2MotMio, (27) 

where the sum over i wil l ag ain be truncated at ^ = 4 in our calculations. 

In earlier publications jiol 20l . 22| . we showed how the computational effort can be reduced 



by expanding the matrix Mqo{E) around the infinite volume energy. Results from these 
calculations agree with those obtained using the full energy dependent matrix as long as the 
shift from the infinite volume energy is smaller than 20%. However, since these shifts are in 
general much larger for volume sizes typical for present day lattice calculations, this earlier 
approach is not applicable in the present work and therefore not followed further. 



III. RESULTS AND DISCUSSION 

Using the formalism laid out in the previous section, we now present the negative energy 
spectrum in finite cubic volumes of varying side lengths. We provide a detailed discussion 
of the results with emphasis on the applicability of the perturbative approaches. For conve- 
nience, the dependence of the energies on the boson mass m is reinstated in the following. 

In this publication, we focus on systems with a positive two-body scattering length a. In 
such systems, a physical diboson state with the energy E^ = — l/(ma^) exists in the infinite 
volume. This energy constitutes the threshold for the break-up of a triboson into a diboson 
and a sing le boson. The volume dependence of the diboson in a finite cubic box is derived 
m |33[ |. The threshold for the break-up into three single bosons is given by E = 0. 



We investigate four states with different energies in the infinite volume: 

• State A: E~ = -10/{ma'^) 

• State B: E^^ = -5/{ma^) 

• State C: E^ = -3/{ma^) 

• State D: E^^ = -1.5/{ma^) 
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In these "universal units" , the energy of the three-body state is given in units of the two-body 
bound state, which is physical for a > and virtual for a < 0. 

For each of these states, the finite volume spectrum has been calculated for several values 
of the box side length L. In applications of the formalism to Lattice QCD, the volume under 
consideration will be smaller than the unnaturally large scattering length. Therefore, we 
concentrate on the region L < a and show results for various box side lengths between L = a 
and L = 0.25a. 

In order to provide an understanding of the universal units employed in this publication, 
we take a look at the nucleonic sector. There, two spin-isospin channels with associated 
scattering lengths ai^^ = —23.7 fm and assi = 5.4 fm are present. Using a "middle ground" 
of 15 fm for the unnaturally large scattering length, we see that the investigated volume range 
contains typical volume sizes of present day Lattice QCD calculations, which are about 2.5 fm 
to 4 fm. The triton energy in units of the deuteron binding energy is Et = —S.SEd, placing 
this state roughly between the states B and C. Note, however, that all results presented in 
this paper are for bosonic systems. 

In order to explicitly verify that our results are renormalized, we use different cutoffs A. 
For all cutoffs considered, the three-body force parameterized by H{A) has been fixed such 
that the same infinite volume energy is obtained. If our results are indeed renormalized, 
the results for the different cutoffs should agree with each other up to an uncertainty of the 
order 1 / (Aa) stemming from the finiteness of the cutoff. In lattice simulations, the scale for 
the cutoff is set by the edge of the first Brillouin zone given by n/b, where b is the lattice 
spacing. The cutoffs employed in this calculation, Aa = 314 and Aa = 225, correspond to 
typical present day lattice spacings of 6 = 0.1 fm and b = 0.15 fm, respectively. 

We investigated five different variants. The Grid approach from Eq. (|5]) is the exact 
result. All other runs employ the Sum approach. The leading order of our perturbative 
approaches is obtained by only taking S-wave contributions into account (cf. Eq. (!26|) ). The 
three remaining runs include the G-wave {i = 4) using the full equation ( !T2|) as well as the 
two perturbative approaches (cf. Eqs. (IT^ . i^S^) described above. 

The volume dependence of the investigated states is depicted in the left panels of Fig. [2] 
for the Grid approach as well as the S-wave and full G-wave calculation in the Sum approach. 
The infinite volume energy of the state is shown as a star on the left axis. The curves for 
different cutoffs agree with each other to within better than 1%, indicating proper renor- 
malization of the results. There is no visible difference between the three depicted variants 
for all investigated volumes. 

The curves show the qualitative behavior known from our previous work for box side 
lengths down to L = 0.5a When decreasing the volume size from infinity, the 



state remains unaffected at first. Going below a certain box size starts to affect the state 
strongly, which is in line with the expected exponential behavior of finite volume corrections 
for bound states. The box size below which the state shows strong deviations from its infinite 
volume energy is tied to the spatial extent of the state in the infinite volume, which can be 
estimated as {—mE^)~^^'^. The size of the state is also correlated with the length scale in 
the exponential behavior of the finite volume corrections. 

To see that this is indeed the case, we summarize the finite volume shifts at L = a calcu- 
lated using the Grid approach for all four states and the two investigated cutoffs Table HIl 
The results for the two cutoffs agree with each other, again showing that the finite volume 
results are renormalized. Comparing the shifts, we notice that state A is shifted by 45%, 
while the more shallow and, therefore, larger state B is already shifted by 137%. For all 
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State 


Efma^ 




A „ OOC; 

iV(2 — Z/0 

E-^(L = a) ma? 


A „ Q1 /I 

E'^{L = a) ma^ 


drol 


A 


-10 


0.32 a 


-14.49 


-14.49 


44.9% 


B 


-5 


0.45 a 


-11.87 


-11.87 


137.4% 


C 


-3 


0.58 a 


-10.63 


-10.63 


254.3% 


D 


-1.5 


0.82 a 


-8.28 


-8.28 


452.0% 



TABLE II. Comparison of infinite volume energies and finite volume energies at L = a as calculated 
by the Grid approach. Shown are results for the two different cutoffs along with the relative 
deviation from the infinite volume energy Jrei = {E^{L = a)/E^) — 1. Also shown is the size 
estimate (— mE"^)"^/^. 

states except the deepest state A the relative deviation from the infinite volume energy 
exceeds 100%. It is even larger than 400% for the shallowest state D. The shift at L = a 
is strongly dependent on the size of the state in the infinite volume, and therefore on the 
infinite volume energy. Thus, a given finite volume affects a shallow state with a large spatial 
extent more than smaller, more deeply bound state, as expected. 

The volume dependence of the three-body energy as shown in the left panels of Fig. [2] 
seems identical for the three plotted approaches. However, there are small differences be- 
tween the different runs. We use the exact result from the Grid approach as a baseline and 
also present the differences A£'3(L) = E^^'^\l) — E^^^^'^\l) in the right panels of Fig. |2J 
Here, E^^™^ can be any of the four runs employing the Sum approach. Namely, these are the 
S-wave only calculation, the full calculation including the G-wave and the two perturbative 
approaches. 

The differences to the respective baseline appear cutoff independent. Some variation 
is visible in the Strictly Perturbative treatment, but it is within the uncertainty from the 
finiteness of the cutoff which is of the order 1 / ( Aa) . 

The full results with and without the G-wave as well as the partially resummed results 
are always larger than the Grid results. The S-wave only result in general deviates from the 
Grid result by more than 0.001/(ma^). The inclusion of the G-wave consistently reduces 
the difference to the full Grid result by about two orders of magnitude. In other words, the 
inclusion of higher partial waves lowers the energy of the investigated state. 

For the very deeply bound state A, the difference between the respective Sum approach 
runs and the Grid approach does only weakly depend on the volume. For volumes with 
side lengths larger than L = a/2, state B shows the same qualitative behavior. Below 
L = a/3, the difference between the Sum approach and the Grid approach becomes larger, 
but is still very small compared to the energy of the state itself. This can be understood 
as higher partial waves with i > 6 having a larger and larger contribution. In addition, 
the absolute size of the G-wave contribution strongly increases, as expected, when going to 
smaller volumes. 

Turning to state C, we observe a behavior similar to that of state B. The strong increase 
in the spread between the Sum approach runs and the Grid approach as well as in the size 
of the G-wave contribution starts at a larger volume, namely around L = 0.6a. This reflects 
the larger spatial extent of the more shallowly bound state C compared to state B. The 
differences of the Sum approach to the Grid results grow by two orders of magnitude within 
the investigated volume range. For the full result, for example, the difference increases from 
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FIG. 2. Left: Variation of the three-body energy with the box side length L for two cutoffs. 
Plotted are the Sum approach with S-wave only and including the G-wave as well as the Grid 
approach. The infinite volume energy is marked as a star on the left axis. The dashed line 
is the two-body energy. Right: Difference between the Grid approach and the Sum approach 
results (S-wave only and including the G-wave via full calculation, partial resummation and strict 
perturbation theory). Note the different volume range in the last row. 
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lO-V(ma^) at L = a to Q.QQA/{ma'^) at L = a/A. For the very shallow state D, we observe 
a strong increase of the differences when going towards smaller volumes. They also grow by 
two orders of magnitude, albeit over a much smaller volume range. The difference of the 
full result at L = a is 2 X 10~^/ (ma^), about twice as large as for state C. It has grown by 
two orders of magnitude to 0.004/(ma^) already at L = 0.43a. 

The size of the higher partial wave contributions is comparable to the uncertainty due 
to the finiteness of the cutoff, which is of order l/(Aa). This should change when including 
higher orders of the EFT. The next-to-leading order corrections are governed by the effective 
range and are of the order r^/a and kr^, where k is a typical momentum. Their inclusion 
will be part of the extension of this framework to the three-nucleon system. 

Overall, the G-wave and higher corrections are indeed small compared to the S-wave 
energies, justifying their perturbative treatment. In particular, the partial resummation 
approach captures the G-wave corrections almost completely. More than 99% of the shift 
from the S-wave result to the full result are accounted for using this approach. 

For sufficiently large box sizes, states C and D show a different volume dependence than 
the more deeply bound states (bottom left panels of Fig. |2]). For both states, we find energies 
above the two-body energy depicted by the dashed line. For state C, this is the case for 
the smallest investigated volume L = a/4 only. The data points for L = 0.3a suggest a 
flattening out of the volume dependence, allowing the two-body energy to catch up when 
further decreasing the box size. This is more clearly seen in the volume dependence of 
state D. At L = 0.6a, the two- and three-body energy in the finite volume are equal. For 
smaller volumes, the three-body state is found above the two-body state. Decreasing the 
volume further increases the energy of the state, until it reaches the three-body breakup 
threshold E = for a volume with side length L = 0.43a. We also observed this behavior 



in the framework developed earlier |20l . |22 



Our interpretation is that, as the two-body energy reaches the three-body energy, the state 
behaves like a boson-diboson scattering state when further shrinking the volume. Indeed, 
the data points above the threshold are compatible with a power-law behavior that is known 
to describe two-body scattering. Note, however, that we are not in the Liischer limit L 3> a 
and therefore do not expect a leading 1/L^ behavior [13]. The presence of three-body effects 
in the volume dependence of such a state will be the subject of further study. It was recently 
shown that, in this topological phase is to be expected |34l |. 



The volume size at which this crossing happens is again tied to the infinite volume energy. 
State D has a smaller infinite volume binding energy, and therefore a larger spatial extent. 
Thus, the box which is not able to hold the complete three-body bound state is larger than 
for the smaller state C. 

Turning back to the differences shown in the right panels of Fig. [21 we note that the 
volume dependence of the higher partial wave contributions does not show any distinctive 
features at the crossover box sizes. 



In the following, we compare the Poisson approach from Refs. [20|-|22| to the approaches 
developed in the present paper. The values for the Poisson and Sum approach in Fig. [3] have 
been obtained by truncating the partial wave expansion of the angular dependence at i = 0. 
We observe a visible spread between the Sum approach and the Grid approach on the one 
hand and the Poisson approach on the other. In the Poisson approach, the sum equation ([2]) 
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FIG. 3. Volume dependence of state B for two cutoffs, Aa = 314 (lines) and Aa = 225 (symbols). 
Shown are the exact result from the Grid approach as well as the S-wave only results from the 



Sum approach and the Poisson approach. The infinite volume energy £"3 
a star on the left axis. The dashed line is the two-body energy. 
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is rewritten as 



Stt 



J2 Z{E-p,q)T{E-q)m 



(Poisson) 1 



(28) 



r 



{E-q):F{q). 



The additional angular dependence introduced by the Fourier transformation leads to a 
recoupling of the partial waves in this approach. We note again that the S-wave component 
of the Sum approach almost completely captures the full result, as discussed before. In the 
Poisson approach, in contrast, the higher partial wave contributions are sizable. 

This is encouraging in view of the extension of this framework to the three-nucleon system, 
which will be the subject of a future publication. In the old framework, we compared the 
volume dependence of the triton 21, 22 1 to results of a Lattice Chiral EFT calculation (23 |. 
The decline of the energy observed in the latter was much stronger. The results described 
above indicate that this discrepancy will be at least diminished when extending the present 
framework to the three-nucleon system. 

We close this section with a discussion of the numerical scaling of the presented ap- 
proaches. The most costly operation in this framework is finding eigenvalues, which is 
equivalent to performing a matrix diagonalization. This is in general an 0{N^) process for 
real non-symmetric matrices of dimension N x N. In both the Grid and the Sum approach, 
the dimensionality of the matrix is dependent on the box side length as it determines the 
size of the discretized momentum space. The left panel of Fig. H] shows the number of vectors 
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FIG. 4. Left: Number of vectors q G ^(Z^) with |^ < A and number of distinct absolute values 
for the two investigated cutoffs. Right: Length dependence of the runtime for three runs and two 
cutoffs. 



q G ^(Z^) with |g| < A, relevant for the Grid approach, as well as the number of distinct 
absolute values of these vectors, relevant for the Sum approach, for the two investigated 
cutoffs Aa = 225 and Aa = 314. As L increases, the vectors in the discrete space ^(Z^) 
lie closer to each other, eventually resembling the continuum. We therefore expect the total 
number of points to be proportional to the volume of a sphere, scaling as (AL)^. This is 
indeed the case for L > a/2. For smaller volumes, the behavior is slightly different because 
the discretized space is still too "coarse" . The number of distinct radii can also be fitted to 
a power law, scaling like (AL)^. 



All results presented in this publication were obtained using a sequential implementation. 
The runtimes range from seconds for the smallest volumes to days for the larger volumes 
in the grid approach. The length dependence of the runtime is depicted in the right panel 
of Fig. m for the two cutoffs and three different runs, namely the Grid approach as well 
as the S-wave only and the full run in the Sum approach. The runtime of the Strictly 
Perturbative approach is almost identical to that of the S-wave run, while the runtime of 
Partial Resummation approach is comparable to that of the Sum approach including the 
G-wave. The scaling of the Sum approach is much better than that of the Grid approach, 
as expected. In the latter, the dimensionality of the matrix that has to be diagonalized in 
an 0{N^) operation scales like (AL)^ compared to (AL)^ for the Sum approach. Therefore, 
we expect the runtime of the Grid approach to scale with three additional powers of AL. 
Indeed, fitting power laws to the runtime behavior yields a scaling of O ((AL)^) for the 
Sum approach and O ((AL)^°) for the Grid approach. The inclusion of the G-wave does not 
change the scaling of the code but merely increases the runtime by about 10%, independent 
of the cutoff. This is a significant improvement in comparison to the Poisson approach 
employed in previous publications, where the inclusion of the G-wave increased the runtime 



by a factor of 4 [20|, l22|. By using parallelized code, one should be able to significantly 
reduce the runtime for larger volumes where the matrix dimensionality exceeds about 2000. 
The volume range where this is the case depends on the approach employed as well as the 
cutoff chosen. 
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IV. SUMMARY AND OUTLOOK 



We studied the volume dependence of three-boson states below the three-body break-up 
threshold inside a cubic box with periodic boundary conditions. A new framework has been 
developed that provides access to the region above the boson-diboson break-up threshold 
and allows for a significantly faster numerical implementation of the higher partial waves. 
We derived an infinite set of coupled equations for the partial waves of the boson-diboson 
amplitude. These equations were solved for several cubic volumes with side lengths ranging 
from L = a to L = a/4. Proper renormalization of the results was explicitly verified. 
We studied the effects from higher partial waves and found them to be comparable to the 
variation due to the use of a finite cutoff. We showed how these contributions can be treated 
perturbatively, either by partial resummation or by the use of perturbation theory for the 
eigenvalue equation. Both approaches yield results that are in good agreement with the full 
calculation, with the former performing slightly better and the latter offering a significant 
reduction of runtime. This is an improvement compared to our previous framework where 
higher partial waves play a more pronounced role and their numerical treatment is more 
tedious. 

The next step is to extend the framework to the three-nucleon system and include higher 
orders of the EFT. The latter is necessary in order to perform precision extrapolations of 
finite volume results. Work in this direction is in progress and will be the subject of a future 
publication. With the new framework, it is also possible to study states that extrapolate to 
a scattering state in the infinite volume. This also allows one to examine, for example, three- 
body effects in Lattice QCD simulations of nucleon-deuteron scattering. In related work, it 
was shown that the compositeness of the deuteron yields modifications of topological nature 
to Liischer's formula j34|. 

We note that an extension of the Liischer formula relating the infinite volume scattering 



phase shifts to the discrete energy levels in a finite volume [10[ is implicitly contained in our 
work. It provides the framework to determine the low-energy constants of pionless EFT in 
the two- and three-body sector from discrete energy levels in a cubic box. After this has 
been done, the infinite-volume scattering observables can be calculated in pionless EFT. 

In summary, our results demonstrate that the finite volume corrections for systems of 
three identical bosons are calculable and under control. The role of higher partial waves was 
assessed and their effect was shown to be of highly perturbative nature. With high statistics 
Lattice QCD simulations of three-baryon systems within reach [I, [l2|, the calculation of the 
structure and reactions of hght nuclei appears feasible in the intermediate future. 
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Appendix A: Derivation of Eq. 1^ 

The partial wave components of the one-boson exchange part of the interaction ker- 
nel Z{E;p,q) are given by the quantity lf\q) defined in Eq. ([7]), 



.c<^o 



a + p-Cq a—p-Cq 



The integration can be carried out by expanding the kubic harmonic in terms of spherical 
harmonics, transforming the integration variables by a rotation R defined by qRcz = Cq 
and using Wigner D-matrices: 



c^om J 4ti \a+p-Cq a-p-CqJ 



CeO m 
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a 



y2y2Cem^Yem{Cq)Ne—Qe . 

-Qi (-) E^^i^(^^) 



Ceo 



pq \pq 



(Al) 



Here, Ni = y is the normalization constant of Y^o- The Wigner D-matrix is given by 

-^om(/^'7'0) = -^Y£m(7,/3), where /3,7 are the Euler angles of the rotation under consid- 
eration, which are by the definition of R just given by Cq. Further, we use a well-known 
integral representation of the Legendre functions of the second kind jssj and, in the final 
step, the invariance of the K^ii under cubic rotations. 
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